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Abstract 

In calculations of the hyperfine splitting in charmonium, the contributions of the disconnected diagrams 
are considered small and are typically ignored. We aim to estimate nonperturbatively the size of the result- 
ing correction, which may eventually be needed in high precision calculations of the charmonium spectrum. 
We study this problem in the quenched and unquenched QCD cases. On dynamical ensembles the discon- 
nected charmonium propagators contain light modes which complicate the extr action of the signal at large 
distances. In the fully quenched case, where there are no such light modes, the interpretation of the signal 
is simplified. We present results from lattices with a « 0.09 fm and a ^ 0.06 fm. 

PACS numbers: 12.38.Gc, 14.40.Pq 
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I. INTRODUCTION 



Lattice calculations of the hyperfine splitting in charmonium usually ignore the contributions 
of the annihilation (disconnected) diagrams to both the vector 7/\[r and the pseudoscalar ric states. 
This simplification leads to an error, and our goal is to determine the actual value of the contri- 
butions. Perturbatively, the contribution of these diagrams in charmonium is expected to be small 
due to the Okubo-Zweig-Iizuka suppression, especially for the vector state [1]. However, non- 
perturbative effects, such as the ^7^(1) anomaly [2] and mixing with glueball and light hadronic 
states, might enhance it enough so that it becomes a nonnegligible fraction of the hyperfine split- 
ting. Previous calculations ||3i |4|] using two-flavor gauge ensembles very roughly estimated the 
contribution to be within ±20 MeV. They both confirm that there are significant difficulties in 
obtaining a signal for the disconnected diagrams due to noise, especially for heavy quarks. 

In our work, the charm quarks are simulated with the clover fermion action with tuned to the 
physical charm quark mass. The disconnected diagrams are calculated stochastically with spin- 
and color-diluted sources. Our calculation improves on the previous ones in a number of ways. 
First, we use larger lattice volumes (28^ x 96, 48^ x 144, and 40-^ x 96) and point-to-point (PTP) 
propagators, which significantly improve our statistics and signal-to-noise ratio. Point-to-point 
propagators reduce the relative standard error over time-slice-to-time-slice (TTT) propagators by 
one to three orders of magnitude. Second, our gauge ensembles have much finer lattice spacings 
a. We work with lattices with a ~ 0.09 fm (fine ensembles) and a ~ 0.06 fm (superfine ensemble). 
Table HI gives the parameters of the ensembles. And finally, we employ the unbiased subtraction 
technique |5] in the stochastic estimators used to determine the disconnected correlators. The 
success of this technique depends on the fast convergence of the hopping parameter expansion of 
the clover Dirac operator used in the subtraction. Considering that is still small for the charm 
quark at these lattice spacings, we use the terms of the expansion only up to third order in K^, which 
reduces the standard deviation of the disconnected correlator by an additional factor of about four. 

In this study we attempt to determine the size of the effects of the disconnected diagrams on 
the mass of the r\c only. Our previous studies [|6|, lZD and our current work show that the effect of 
the charm annihilation on the vector state are much smaller than 1 MeV; thus we ignore it here and 
equate the hyperfine effects with the effects in the pseudoscalar only. Our calculations are done on 
two fully quenched and one dynamical ensemble with two light degenerate quarks and one strange 
quark (2-1-1 dynamical flavors) in the asqtad formulation ||8|]. In the fully quenched case, the 
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Ensemble a [fm] m/ / Volume # config. 

QF « 0.085 ••• 28^x96 0.120,0.127 410 

QSF w 0.063 ••• 48^x 144 0.125,0.130 415 
DF « 0.086 0.0031/0.031 40^ X 96 0.125, 0.127 766 

TABLE I: Run parameters of the quenched fine (QF), quenched superfine (QSF) and dynamical fine (DF) 
ensembles are shown. The bold values of are the ones obtained by tuning the rj^ mass and are used in 

lies 1 3, 



this study. The nonbold values are from our previous studies 



70 and are listed for comparison. 



disconnected Tj^ correlator can have at most additional contributions from the C/a(1) anomaly and 
close-lying glueball states. In the 2+1 flavor dynamical case, the disconnected correlator can also 
couple to light hadronic states, which complicates the task at hand significantly. In both the 2+1 
flavor dynamical and fully quenched cases we ignore contributions to the disconnected correlators 
from sea charm quark loops . To the extent that the disconnected contribution is small (first order), 
the sea charm quark effects are second order and so negligible at our level of precision. Our result 
that the contribution is, indeed, small makes the calculation self consistent. 

This paper is organized as follows. Section II outlines the analytic framework in which we 
interpret our lattice data. Section III discusses some general properties of the disconnected propa- 
gators which follow from our analyses. In Sec. IV we give our fitting method for the disconnected 
propagator. Section V is dedicated to the specifics of the tuning of the charm quark mass. The 
final section, VI, contains our results and conclusions. 

II. GENERAL ANALYTIC FRAMEWORK 

In this section we derive the shift of the mass of a flavor singlet state due to the contribution of 
the disconnected diagrams to its full propagator. Figure[I] shows the diagrammatic expansion of the 
full propagator, where with continuous lines we represent the charm quark propagators. The first 
term in this expansion is the connected piece and the rest are disconnected propagators containing 
charm quark loops. We denote the momentum- space connected propagator of a (pseudo)scalar 
meson as 

C(/) = , (1) 
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where A is a constant, and is the "connected" mass of the meson. (The vector meson propagator 
has the same form as Eq ([1]), if we neglect the spin degrees of freedom). Then the full propagator 
is the infinite sum: 



= 2^ 2 + 2^ 2 ^ P 2^ 2 + 2^ 2 ^^ 2^ 2 ^^ 2^ 2 + ' ' ' ' 



where the first term is the connected piece C{p^) and the rest are terms which consist of dis- 
connected quark loops of the same flavor as the meson's constituent quarks (see Fig. [U for the 
diagrammatic representation of the three explicitly given terms). The function ^{p^) effectively 
describes all possible interactions between the quark loops in the disconnected pieces and the 
gauge fields, quarks of other flavors or effects such as the t/A(l) anomaly, if relevant for the spe- 
cific meson state. The disconnected propagator D(p^) is naturally the sum of all terms in Eq. 
except the first one. After we sum the geometric progression in Eq. (O, we obtain: 

where my^ is the "full" meson mass we could calculate if we knew all terms that contribute to the 
full propagator. Thus, the difference between the mass that is usually computed from only the 
connected propagator and the actual mass is approximately ^ 

om = nic — mf ^ . (4) 

In the last expression for simplicity we replaced the function ^{p^) with the value of its largest 
contribution in Eq. (O at the pole p^ = —m^. Hence the sign of 5m depends on the sign of m^). 
The mass shift dm due to the disconnected quark loops can be treated as a perturbation, in which 
case, to first order, both the connected and disconnected contributions can be computed without 
dynamical sea quarks of the constituent's flavor (heavy-quark-quenched case). In this case, only 
the first two terms in Eq. ^ survive and the disconnected propagator is reduced to the second term 
only (shown diagrammatically in Fig.[Tl middle). 



' The momentum dependence in X{p^) leads to a first-order (in A,) renormalization of the connected pole residue, so 
there is a second order correction in our result that we safely ignore. 
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FIG. 1: Connected and disconnected diagrams contributing to the full propagator on lattices unquenched 
with respect to the charm quark are shown. 

III. PROPERTffiS OF THE DISCONNECTED PROPAGATOR 

The asymptotic behavior at large times t of the full charmonium propagator, F{t), is 

F{t)=C{t)+D{t) = Y^{Q.\0\n){n\0\Q)e-^''^ (£2|C>|0) V'^"' , (5) 

where the sum is over all eigenstates \n) of the Hamiltonian with corresponding energy eigenvalues 
En and is the vacuum state. In the last part of the above expression Eq is the mass of the lightest 
state contributing to F{t). The operator O is defined to be Hermitian, in which case F{t) > for 
all t. This is also true if we consider the PTP propagator F{r) instead, where r is the Euclidian 
distance. The matrix defining the spin structure in the operator O is F = iy^^iyi for the T|c and 
states, respectively, in terms of Hermitian 75 and Yj. At large distances r, the lightest possible 
modes that couple to the operator O should dominate m F{r). The origin of these can be light 
glueballs and, in the dynamical case, the propagation of hadronic modes consisting of quarks 
lighter than the charm quark. Since F(r) is nonnegative for all r, it follows that, when it dominates, 
D{r) should also be nonnegative in the large distance limit. The sign of D{r = 0), with the above 
hermiticity condition on O, is strictly negative for the pseudoscalar (and positive for the vector). It 
follows that in the dynamical case, where D{r) is dominant at large distances, D{r) changes sign 
for the pseudoscalar. In the quenched case this sign flip occurs if there are glueballs lighter than 
the charmonium state studied. On the lattice, the TTT disconnected propagator is calculated as: 

D{t) = cr{L{OyL{t)), where L(?) = Tr(rM-i), (6) 

and M is the charm quark matrix. The trace in L{t) is over the Dirac, color and space indices. For 
the vector we have T = iji, cr = 1 and for the pseudoscalar T = ijs, cr = —I. On the other hand, 
the PTP discoimected propagator is obtained in the following maimer: 

r=\x-y\ 
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FIG. 2: The comparison of D{t) and D(r) for r\c for two different ranges of t and r is shown. The results 
are from the calculation on the dynamical ensemble at = 0.127. The data for D(r) was averaged in small 
bins in r for r > 5. 



where x and y are lattice coordinates, the sum is over all pairs of lattice points at Euclidean distance 
r = \x — y\, Nr is, the number of these pairs and there is a trace only over spins and colors in L{x) 
(but not over space). From the previous studies M>\^ it is known that the D{t) signal disappears 
very quickly around t = 2 — 3. We work with the PTP disconnected propagator instead, since 
this way we benefit from both the additional data at noninteger distances and the much improved 
statistics. The correlator D{r) has from one to three orders of magnitude smaller relative errors 
than the TTT disconnected propagator in the region where we have a signal. Figure|2]illustrates this 
statement by comparing D{r) and D{t) for r[c for two different ranges of r and t. Both propagators 
are calculated with Kf = 0. 127 on the dynamical lattices from Table I. In the right panel of Fig. |2l 
the comparison is done on a shorter range in order to emphasize the fact that we do have a clear 
signal for D{r) in the range where the D{t) signal is completely obscured by the noise. The result 
that the D{r) signal is so much better than the one for D{t) can be explained by the fact that in 
the calculation of D{t) there are a great number of contributions from points, which, although not 
far from each other in the t direction, are far in the 4d Euclidean space. For the disconnected 
correlator, the noise increases strongly with the distance and such points contribute nothing to 
D{t) but noise. This problem is avoided when working with D{r) instead. We also note that as 
predicted in the previous paragraphs, both D{r) and D{t) undergo a sign flip for the r|c state. 
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IV. FITTING THE DISCONNECTED PROPAGATORS 



To determine 5m for the r\c we use Eq. (U), which means that we have to obtain from 
our data for the PTP disconnected propagator. In order to fit our data for D(r), we need a fitting 
model which satisfies the requirement that the charmonium disconnected propagator is treated as a 
composite object, which has contributions not only from the studied charmonium ground state, but 
also possible effects from excited charmonium states, states lighter than the charmonioum ground 
state, and possibly the UA{i) anomaly. We also have to take into account that our data exhibits 
rotational symmetry violations at short distances, due to the finite lattice spacing. 

To define such a fitting model it is easiest to start from the momentum-space description of the 
disconnected propagator. A simplified form which describes its behavior in momentum space is 

Dip') = Hp') (^^2 + 1 ^^^] ' (8) 

where we have included in the quark loops one ground state, characterized by mass and A'^ 
excited states with masses m" (the index n= 1, . . . , A''). Here we also make the assumption that the 
interactions for all states are described by the same function 'k{p^). In the fully quenched case, we 
model the function as 

X{p^) = U + , (9) 

where U stands for possible effects of the L'a(I) anomaly, so it is negative, and f / {p^ + m^) is 
an effective light glueball term with lUg being the glueball mass. We assume that U, f. A, and 
^n=i...N cj^ange little for a wide range of momenta and will approximate them with constants in 
our model. In the 2+1 flavor dynamical case, the expression for could be more complicated. 
For example, we need to take into account the existence of both light glueballs and light hadronic 
modes in order to describe our data: 

X{p^)=U + ^^ + ^^ + ---. (10) 

In the above / is a constant and m/ is the mass of one of the light hadronic modes. In practice we 
keep only one light hadronic mode in the above with an effective mass m; that (we hope) describes 
well the long distance behavior of the PTP propagator. 

We want to limit the number of free parameters in our model to as few as possible, since 
although our data is of much higher quality than in other studies, it is still difficult to resolve all of 
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the parameters in Eq. ([8]) from the disconnected propagator data. The masses m^- and m^^^- for 
example, can be determined from fits to the TTT connected charmonium propagator C(t) , and then 
be used as constants in our fits. We obtain the ground state mass rric with a very small fitting error, 
but the excited states masses m"^^ -^ are less well known. We can also determine with varying 
degrees of precision the constants A and A"=^■^^^ from fits to C{t), since they are proportional to 
the amplitudes of the ground and the excited states, respectively. We relate the PTP amplitude A to 
the corresponding TTT amplitude as follows: We Fourier transform the ground C{p^) propagator 
and then integrate over space: 

The right-hand side in the above is the TTT propagator, which implies a relation between the 
amplitude of the ground state in C{t) = AfC^"^^^ H , denoted by At, and the factor A: 

A = InicAt . (12) 

Another parameter that we can fix in our model using prior knowledge is the glueball mass mg. 
We use the results for the lightest ^ glueball from Ref. [|9|], namely we set mg = 2563 MeV, 
which is the value extrapolated to the continuum limit. We use this value for all of our lattice 
spacings, since in Ref. [|9|] it was found that the glueball mass does not vary much at fine lattice 
spacings and is compatible with the continuum extrapolated result. On the dynamical ensembles 
we have to take into account also the contribution of the light hadronic modes, and preferably we 
want also to set the mass m/ to an appropriate constant. In our previous work [la LZD we found 
that the long distance behavior of the PTP pseudoscalar propagator on the dynamical ensemble 
can be fitted well with a light state of mass ami ~ 0.42. This is very close to the physical mass 
of the Tj' of 958 MeV; thus we fix m/ to the mass of the T]'. Although there are states lighter than 
the Tj' contributing as well (such as the T], multipion states etc.), this approximation is probably 
satisfactory, because the modes lighter than the Tj' would mix even less with the heavy Tjc state, 
and we cannot distinguish their signal at our level of statistics. Thus the only free parameters in 
our model remain U , f, and / (the last one is present only in the dynamical case). 

The summary of our fitting strategy in the simpler fully quenched case is as follows: 

• On a given lattice ensemble we calculate the TTT connected propagator of the r\c state. 
From fits to it with the asymptotic form C{t) = AtC^'"'^' +I^^^jA"e^™"' we determine m^. 



m1^^-^. At, and A;'=i Using Eq. ^ with the substitution mj"^ ^ sj 2(cosh(mJ."^) - 1) 
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in order to take into account lattice discretization effects, we obtain A and A"~^ We use 
the central values of all of the above parameters as constants in our model function Eq. ([8]). 

• We also fix the parameter nig in Eq. using prior knowledge. 

• In Eq. dH]) we replace with E(2(l — cos(j9,)) and all the masses with the appropriate 
expression as done in the first item above for m["\ to account for the lattice discretization 
effects. Equation ^ can be rewritten in a form linear in the two parameters U and /, which 
is more convenient for fitting purposes: 

D(/) = [/ri(/)+/r2(/). (13) 

Next, we do a discrete Fourier transformation of Ti2{p^) on a lattice of an appropriate size 
and obtain the functions Ti^^r) at discrete values of r. We tabulate Ti^^r) at each distance 
r < 15. This range of r is sufficient, since our signal is too noisy for r > 15. Thus, using the 
linear model 

Dfit(r) = [/ri(r)+/r2(r), (14) 

we fit our data for the PTP disconnected charmonium propagator D{r) in order to extract U 
and /. 

• With the fit values of U and / at hand, we determine X{—m^.) and 5m from Eqs. (H]) and 

The fitting procedure in the 2+1 flavor dynamical case is quite similar. In addition there we fix 
appropriately the light mass m/ and use a fitting form with three tabulated terms: 

Dmir) = UT,{r) + fT2{r) + l^r) . (15) 

After we extract the parameters U, f, and /, we solve again for 5m with the appropriate X{—m^) 
fromEq. 

The error due to the assumption that the participating masses and the amplitudes A and A"^^- 
in our model are constants is discussed in Sec. VI. The success of our fitting model depends on 
how well it approximates the interaction dynamics on the lattice and on the quality of our data. 
An essential part of the construction of the fitting model for the PTP disconnected propagator 
turned out to be the number of excited states A'^ which we have to include. The excited states 
have larger contributions to the disconnected propagator at a given distance than they have in the 
case of the connected propagator. This means that a good knowledge of the spectrum of excited 
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states is required in order to fit the disconnected propagator. In the disconnected diagram with 
net propagation distance r the connected r[c propagates on average less than half the distance 
before annihilating and again less than half after reappearing. This means that a reliable fit to the 
connected propagator yields masses and amplitudes which we can use in fits to the disconnected 
propagator at least at twice the distance, where the signal may be too noisy. Thus, we are limited 
to fits of the connected propagator with ^min = 2 or 3 and use the extracted parameters in the 
disconnected fits from r^in = 4 or 6. 

To obtain the charmonium spectrum from the fits to the connected propagator, we employ a 
fitting model which forces the splittings between the states to be positive, essentially creating a 
"tower of states". The priors for the logarithms of all the splittings are the same {i.e., we assume 
the states are equidistant as a first approximation) and so are their widths. We used a set of different 
values for the splitting priors and their widths to check the stability of the resulting spectrum. We 
found that the extracted masses were stable for a relatively wide range of priors. Our best fits have 
the splitting priors in the range of 570-770 MeV. Still, this approach is not intended to provide a 
reliable determination of the masses of excited quarkonium. Rather, we use this heuristic model 
to understand the role played by excited states in our result. 

In our fitting model we also require (through the use of priors) that the amplitudes of the excited 
states are no larger than the ground state amplitude. Without this restriction, we noticed that the 
amplitudes of the excited states grow noticeably larger than that of the ground state (where we 
define this to mean difference larger than 1.5o). The assumption that this should not happen 
stems from general considerations: The wave function at the origin is smaller for an excited state 
than the ground state, simply because the excited state wave function spreads out more. This is 
a characteristic of nonrelativistic potential models. There is also some support from experiment: 
The decay constant of ^]f' is smaller than that of the J hi [PDG values: 279(8) MeV 411(7) 
MeV, respectively]. Lattice studies of the light-light lIlOll and heavy-heavy [11] meson sectors 
also confirm this expectation. But other lattice calculations, such as in Ref. [|l2l] for the heavy- 
light meson case, show the first excited state decay constant growing larger than the ground state 
one. A possible explanation for this discrepancy in the last study is that the contributions of the 
neglected higher excited states became lumped into the amplitude of the first excited state. 

In our fitting model, the requirement that the amplitudes of the excited states do not grow much 
larger than the ground state amplitude, we hope, prevents the "clumping" of states with similar 
masses into an effective state with a large effective amplitude. However, since we use only a 
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FIG. 3: (Left panel) Relative positions of the charmonium states with respect to the 0^+ glueball mass with 
nig = 2563 MeV on the quenched superfine ensemble, for different kc, are shown. (Right panel) Same for 
the dynamical ensemble. 

Gaussian prior to constrain the logarithm of the excited state amplitudes, it still often happens, 
depending on the number of states included, that the highest state in the fit and sometimes the 
next highest end up with large amplitudes. We interpret this outcome as a clumping of multiple 
unresolved states of similar mass. How we count them affects our result for 5m. To compensate 
for this effect we represent the resulting contribution as a sum of states of similar mass: 

^ A' 



M 

E 

yt=l 



(16) 



where M ^ A^ /A' = Af' /At (the ratio of the highest excited state amplitude to the ground state 
one). Each of the states in the sum above should contribute a/A^/M/ (p^ + (m^)^) in Eq. ([8]). 
Since there are M of them, the contribution of the effective A^-th state in Eq. ^ is modified to 



VA^M I (/)^ + (m^)^). In effect, we multiply the amplitude by the number M of its contributing 
states before using it in Eq. ([8]) in order to correctly account for the possibility that the A^-th state 
is an effective one. Note that increasing the multiplicity in this way always decreases 5m. Thus in 
our analysis of systematic effects, we explore the sensitivity of our result to this assumption. 



V. TUNING THE CHARM QUARK MASS 



As we already stated in the Introduction, we use clover fermions to generate the connected and 
disconnected charmonium propagators. This means that we have to tune the hopping parameter, 
Kc, to correspond to the physical charm quark mass. In our preliminary work [laLZl] we tuned Kc 
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using the kinetic mass of Ds i.e., we used the Fermilab interpretation of the clover fermions. For 
the fully quenched ensembles this tuning was rather approximate. In this work we adopt a different 
approach, namely, we tune by matching the rest mass of the Tjc to its physical mass instead with 
an accuracy of several percent. Our current approach is more appropriate for the purposes of 
determining 5m, since it positions the charmonium r[c state correctly with respect to the lightest 
glueball with which it might mix. This is important, because, depending on whether the mass of 
the charmonium state is heavier or lighter than the lightest glueball, 5m might change in absolute 
value or even undergo a sign flip. Figure [3] illustrates this statement by showing the masses of 
the r|c, T]'. and the lightest ^ glueball on the superfine quenched (left) and the dynamical fine 
ensembles (right), for values of Kf obtained by our previous and current tuning methods. In the 
quenched superfine case, for example, using the value K^. = 0. 130 from the kinetic mass tuning of 
Ds, gives a rest mass of the Tjc lighter than its physical value and lighter than the lightest glueball. 
This implies that their mixing will "push" the V[c mass to lower values on the lattice. However, 
if the Tjc rest mass assumes its correct physical value (achieved at Kf = 0.125), the effect of the 
glueball mixing would be exactly the opposite. We conclude that although the kinetic mass 
tuning is the correct method in cases when we want to determine various mass spectrum splittings 
in the charmonium system, for our study the appropriate method is to tune the charm quark mass 
using the rest mass of the charmonium state we are interested in. (An alternative method which 
could render both the charmonium splittings and the rest masses correct, is to have different values 
for the spatial and temporal hopping parameters, a strategy which we do not employ here.) 



VI. RESULTS AND CONCLUSIONS 

In this section we present our results for the mass shift 5m due to the contribution of the discon- 
nected diagrams for the Tjf on all of the ensembles from Table HI We use 72 Z2 random sources per 
lattice with spin and color dilution to compute the disconnected propagators on all of our ensem- 
bles (which means there are 72 x 12 quark matrix inversions per lattice performed). To explore 
the systematic effects which arise in the determination of 5m due to our incomplete knowledge of 
the charmonium spectrum, we studied in more detail the data from the quenched fine ensemble. 
We fitted the connected propagator using 4, 5, 6, and 7 states and used the extracted masses and 
amplitudes to fit the disconnected propagator Dr|^.(r), as described in Sec. IV, in each case. 

Table |ll] shows the results from the fits to the connected propagator in its upper and middle 
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#states ffJri^ ml m^ ml m^ m^ m^ X^/DOF 

4 1.3708(2) 1.67(3) 1.98(6) 2.59(5) 1.6 

5 1.3708(2) 1.68(3) 1.98(7) 2.56(5) 2.81(13) 1.6 

6 1.3708(2) 1.67(3) 1.94(6) 2.16(9) 2.42(10) 2.67(11) ••• 1.4 

7 1.3709(2) 1.67(3) 1.93(7) 2.13(8) 2.34(9) 2.54(9) 2.74(9) 1.3 



#states At A} Aj A^ Aj Af Af 

4 0.589(3) 0.58(19) 1.63(34) 5.64(31) 

5 0.589(3) 0.59(19) 1.58(35) 4.81(45) 0.91(36) 

6 0.589(3) 0.56(18) 1.06(33) 0.93(35) 0.86(34) 4.55(50) ••• 

7 0.588(3) 0.56(17) 0.90(31) 0.95(35) 0.95(37) 0.97(38) 3.74(60) 



#states 5m [MeV] 5m™"' [MeV] 

4 -3.61(24) -1.59(15) 

5 -3.31(29) -1.62(15) 

6 -2.74(24) -1.75(20) 

7 -2.45(22) -1.88(18) 

TABLE n: Masses (top) and amplitudes (middle) in lattice units extracted from fits to the connected prop- 
agator on the QF ensemble with different numbers of states (4, 5, 6, and 7) are shown. The lowest part of 
the table shows the mass shift 8m, calculated using the results from the upper two parts of the table. Also 

shown is the "corrected" mass shift 5m^"'^, which is obtained using the systematic error estimation method 
described at the end of Sec. IV. 

subtables. All of the fits were performed on the same time range (t = 2 — 45). We see that the 
X^/DOF improves with adding more states to the fitting model and the amplitudes of the excited 
states become smaller at the same time. The extracted masses of the ground and the two lowest 
excited states appear to be quite independent of the number of states included in the fit. The third 
and higher excited states on the other hand do depend on the number of states. This is not very 
surprising; these states are much more difficult to extract and are likely to be effective states. Their 
amplitudes also are more likely to grow large, lending support to this interpretation. 
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Ensemble m-n ml m? ml mf. ml mf ml mf. 

•[c LLC L L L L L 



QF 1.3708(2) 1.67(3) 1.94(7) 2.16(9) 2.42(10) 2.67(11) 

QSF 0.9734(2) 1.22(2) 1.41(7) 1.63(9) 1.87(10) 2.10(10) 

DF 1.2749(4) 1.58(3) 1.88(7) 2.08(9) 2.26(8) 2.43(8) 2.62(10) 2.83(11) 3.05(12) 





Ensemble 


4 4' 42 43 44 45 46 4? 48 

iv^ A t T t t t t t 


QF 

QSF 

DF 


0.589(3) 0.56(18) 1.06(33) 0.93(35) 0.86(34) 4.55(50) 

0.270(2) 0.32(10) 0.39(13) 0.33(13) 0.33(13) 4.44(24) 

0.822(7) 0.84(21) 1.00(37) 1.30(47) 1.40(49) 1.30(47) 1.17(43) 1.07(40) 0.98(38) 



TABLE III: Parameters extracted from fits to the TTT connected r\c propagator C{t) for each ensemble are 
shown. Masses and amplitudes are in lattice units. The y} /T>0¥ for the fits is 1.4, 2.2, and 1.4 for the QF, 
QSF, and DF ensembles, respectively. The ?niin for these fits is 2, 3, and 2. The central values of the masses 
and amplitudes are used as constants in the fitting model for the disconnected r\c propagator. 

The effect of the number of states in each fit on 5m is shown in the lower part of Table HIl We 
also show there the corrected value of the mass shift 5m'^°'"'^, which is obtained by modifying the 
amplitudes that are significantly larger than the ground state one {i.e., the amplitudes of the second 
and third excited states for the 4- and 5-state fits and the highest excited state one in the case of 
6- and 7-state fits), in the manner described at the end of Sec. IV. The difference between 6m and 
g^corr giyes some idea of the systematic error. This systematic error grows when there are fewer 
states in the fitting model, possibly because the charmonium spectrum is not well represented by 
it over the fitting range. This effect is also signaled by a growing x^/DOF in these cases. It is 
encouraging that 5m'^°'"'^ is very consistent between the different fits. 

Among these results for 5m, the one we consider best is the one obtained using the 6-state fit 
to the connected propagator, for two reasons. First, this is the fit with the fewest excited states to 
achieve x^/DOF below 1 .5, a value we consider to be on the boundary between "good" and "bad" 
fits. Second, the excited state spectrum is close to the picture where all amplitudes but that of the 
highest included state are not much larger than the ground state amplitude. We expect the highest 
excited state to reflect the fact that we work with finite number of states in the fitting model, and 
thus likely to have a large effective amplitude. 

We apply the same criteria when we repeat the whole calculation on the rest of the ensembles 
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in this study. We adopt the results for 5m from the best of our fits as the final answer and report 
the difference between 5m and 5m'^"'^ as an asymmetric systematic error. Results for all three en- 
sembles are summarized in Table Hill The number of states in the preferred fit is chosen according 
to the criteria described above, resulting in six states for the quenched ensembles and nine for the 
dynamical one. For the quenched superfine ensemble however we did not obtain a good x^/DOF 
for ?min = 2 or 3, even when we included more than six states. 

To complete our results for each ensemble, we fit our D-q^{r) data, as described in Sec. IV, 
with the values for the model parameters taken from Table Hill respectively. The glueball mass in 
our model is the constant nig = 2563 MeV, and in the dynamical case we also fix m/ = 958 MeV. 
Figures |4] (both panels) and |5] (left panel) show our results for the r\c disconnected propagator 
for each of the studied ensembles and the best fits to the data. In Table |IV] we give the fitting 
parameters, ranges, per degree of freedom and our final result for the mass shift 5m for each 
ensemble. 

All of our results consistently show 5m < (meaning an increase on the r\c mass). It is no- 
table that this is opposite the prediction of perturbation theory |l3|. Even without a quantative 
calculation of the mass shift 5m one might expect that it will be negative. First, according to the 
mixing models lll4ll . if the Ua{1) anomaly has an effect on the Ti^, it should make its mass larger. 
In our model this is reflected in the fact that we obtain U <0. Second, the light glueball happens 
to be lighter than the Tjc, and when they mix, similarly, the mass of the r\c is pushed up. This is 
reflected in our finding that f/{—m^ + m^) < 0. The same is valid for the effect of light hadronic 
modes on m/ in the dynamical case. The mass shift itself is similar for the two quenched en- 
sembles: 5m = —2.74(24) and —2.18(47) MeV for the fine and the superfine one, respectively, 
where the errors are statistical only. We estimate the systematic effects stemming from our lim- 
ited knowledge of the charmonium excited states as described in Sec. IV, to be around 1 MeV, 
applied in the direction of decreasing the absolute value. This estimate is based on the difference 
between the above values for 5m and their corresponding 5m'^°'"'^, the latter being —1.75(20) MeV 
and —1.15(28) MeV, respectively. The quenching of the light quarks is, of course, another source 
of systematic error which might not be negligible. 

The dynamical ensemble yields a larger value for the r\c mass shift: 5m = —8.52(24) MeV (the 
error is statistical only). We did not determine the systematic effects of the excited states (through 
g^corr-) ^YyIs case, because the amplitudes of the excited states in the connected propagator did 
not grow much larger than the ground state one (i.e, they differ less than 1.5a). This is not too 
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surprising considering that we already use 8 excited states in our fit to the connected propagator in 
order to achieve a good x^/DOF. The discrepancy in the determined values of 5m in the quenched 
and the dynamical cases is most likely due to much larger systematic errors (other than the excited 
states contribution) in the latter case. For example, our simplified model does not account for the 
complications due to the mixing with the open charm threshold. Above this threshold there are 
numerous close-lying open-charm states. These discretized (in finite volume) continuum states do 
not fit the tower-of-states model: At our lattice size their level spacing is approximately ten times 
smaller than the typical level spacing in our tower-of-states model, and they have a degeneracy that 
grows rapidly with S-wave phase space. Further, their amplitudes are likely to be much smaller 
than the amplitudes of the "bound" states. They would appear as clumped, effective states in 
the tower-of-states model. To the extent they are important in our analysis, correcting for their 
clumping would tend to reduce the absolute value of 6m. 

In Table HV] we also present separately the contribution of the anomaly 6mf/ to the total mass 
shift for each of the ensembles. The values of dmy are roughly the same as the one predicted in 
Ref. ImI] using mixing models. Our numerical results show that the contributions in MeV of the 
[/^(l) anomaly is about half to 2/3 of the final value. The effect of mixing of the r\c with light 
hadronic modes in the dynamical case is much smaller than 1 MeV and is practically negligible. 
This is not unexpected considering the large mass difference between them. Figure [5] (right panel) 
shows the values of 5m for each ensemble. 

In addition to the statistical error in these values, another source of uncertainty in 5m is the 
tuning. We have negligible Kf-tuning errors for the superfine quenched (and the dynamical) case 
where the tuning is done to 1% — 2%, but the r\c mass in the case of the quenched fine ensemble 
is about 7% heavier than the physical one. This leads to an asymmetric correction to the 5m in 
the case of the quenched fine ensemble by about 1 MeV in the direction of increasing its absolute 
value. We obtained this correction by assuming that the parameters U and / change negligibly for 
small mass fluctuations and there is only an explicit dependence on nic in Eq. (|4]). Then we equate 
the systematic error with the difference in 5m when we use the physical and the measured value 
of the mass of the r\c. 

Finally, the last source of uncertainty in 5m originates from the assumption that the masses and 
amplitudes in the fitting model Eq. ([8]) are constant, when in reality we know them up to their 
statistical errors. We estimate the effect of this assumption by varying these masses and the ampli- 
tudes within their statistical errors when fitting the disconnected propagator and recalculating 5m. 
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We found that the resuhing error for the quenched ensembles is within the statistical uncertainty 
and we neglect it in the final error budget for these ensembles. For the dynamical ensemble this 
error turned out to be 2 MeV. In this case we add it to the statistical error for this ensemble in 
Fig. |5] (right panel). 

We mentioned in the Introduction that we equate the effects of the disconnected diagrams on 
the hyperfine splitting in charmonium with the mass shift they induce in the r\c only. In other 
words, we ignore the possible mass shift they cause in the J/\\f. We base this approximation on 
our attempt to estimate the effect on the vector using a fitting procedure similar to the one we used 
for the r\c- Our data for the vector is more noisy and the signal in the PTP propagator dies out at 
shorter distances than in the pseudoscalar case, due to the larger mass of the We found the 
effects of the disconnected diagrams for the vector is much smaller than 1 MeV and thus, they are 
within the statistical error of 5m for r\c- 

In conclusion, based on our results for the mass shift in r\c in the quenched case, the charmo- 
nium hyperfine splitting is decreased by 1 - 4 MeV when we take into account the disconnected 
diagrams. This range is represented visually by the band between the slashed blue lines in Fig. [5] 
(right panel). In this final range for 5m we ignore the dynamical result on the basis of its much 
larger and much less reliably estimated systematic effects, which require further study. 
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FIG. 5: (Left panel) The rj^ disconnected propagator vs the EucHdean distance r in lattice units for the 
dynamical fine ensemble at = 0.125 is shown. (Right panel) The mass shift 5m for the QF, QSF, and DF 
ensemble. The black error bars denote the statistical errors for the QF and QSF results. For the DF ensemble, 
the black error bar also includes the error due to specific assumptions in the fitting model as described in the 
text. The asymmetric red error bar originates from the systematic error due to our incomplete knowledge 
of the excited charmonium spectrum. The asymmetric green error bar is due to the mistuning in the QF 
case. The band between the two slashed blue lines encompasses the most likely range of 8m based on the 
quenched data. 
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